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§1. Introduction 

One of most important unresolved issues concerning 
disordered spin systems is whether the low-temperature 
phase of spin glass systems is described by_the mean-field 
pictured or rather by the droplet picture!! For three di- 
mensions, the issue has not yet been settled in spite of 
a lot of arguments and numerical works. W On the other 
hand, most researchers believe that the spin glass model 
in two dimensions with a continuous bond distribution is 
disordered (paramagnetic) at any finite temperature and 
is critical right at the zero temperature. In this case, the 
critical phenomena near the zero temperature have been 
understood mainly by the droplet picturecP and other 
similar arguments such as the domain wall rciiormaliza- 
tion group argument .lP The droplet picture£§ can also 
explain fairly well the experimental results on the two- 
dimensional spin glass systems.© According to the stan- 
dard droplet picture, the characteristic length scale at 
each temperature, i.e., the correlation length, is directly 
comparable to the typical size of an activated droplet. 
Then, it is natural to expect that the various length 
scales associated with the same energy (temperature) 
scale should be essentially the same except for some con- 
stant prefactor regardless of difference in the definitions. 
For example, the correlation length at temperature T 
should be the same as the typical droplet size whose ex- 
citation energy is fcflT. This yields the scaling relation 

-Od = Vt 

where y t = l/v is the thermal scaling exponent, and 9o 
is the droplet exponent that characterizes the size depen- 
dence of the droplet excitation energy. Another example 
is the stiffness exponent 9s that relates the size of the 
whole system to the excitation energy of the domain wall 
which is induced by some appropriately chosen boundary 
condition. The excitation energy is called domain wall 
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energy or "stiffness" and hence the name of the exponent. 
One, then, might expect that yt, —9d and — 9s are equal 
to each other. _In fact, in the domain wall rcnormaliza- 
tion argumentJj'LP it was concluded that y t — —9s- How- 
ever, we see in what follows that a number of numerical 
calculations, such as the direct calculation of the magne- 
tization as a function of magnetic fieldEr and tha-recent 
direct measurement of droplet excitation energies,™ indi- 
cate that 9s is truly different from the other two, whereas 
the estimates of yt and 9d are at least close to each other. 

In this paper, we concentrate our attention on the 
Edwards- Anderson (EA) models described by the follow- 
ing Hamiltonian, 



H 



7 j JijSiSj, 



where Jy is a quenched random variable with some given 
bond distribution. It is a general belief that models de- 
scribed by this Hamiltonian can be classified into two 
classes according to the type of the distribution of J^. 
One is the class of continuous distributions. The Gaus- 
sian bond distribution is the most often studied example 
in this category. Models in this category do not have 
non-trivial degeneracy in the ground states. The other 
class is that of discontinuous distributions to which the 
±J model belongs. Models in this category in general 
have non-trivial degeneracy and the entropy per spin ap- 
proaches to a finite value in the zero-temperature limit. 
In other words, the model is not completely "frozen" 
even at zero-temperature, making the two-point corre- 
lation function decaying to zero. As a result, they are 
characterized by non-zero values of the exponent r/. 

In section 2, we consider the model with continuous 
bond distribution. We review some of previous numerical 
calculations for models in this category focusing on the 
three exponents, &s, y t and 9d- We also present some 
new results on 9s- In section 3, we review briefly recent 
results on the discontinuous distributions. 

§2. Models with Continuous Bond Distributions 



1 



2.1 Domain-wall renormalization- group argument and 
stiffness exponent 
It was argued© that the excitation energy E\y(L) of a 



at zero-temperature was estimated through an extrapo- 
lation of the domain wall free energy at finite tempera- 
tures to zero temperature. The resulting estimate of the 
stiffness exponent was 




Fig. 1. The domain walls in 16 randomly chosen samples with 
L = 48. 



domain wall across a system of size L can be regarded 
as the renormalized coupling constant after the system 
is renormalized up to scale L. It follows that if we define 
the scaling exponent 9s by 

E W (L) oc L 9s , 

9s is expected to equal iu- Also, in the standard droplet 
theory for spin glasses,© E\y(L) is considered to be of 
the same order as the excitation energy of a droplet of 
scale L. If these interpretations are correct, we can es- 
timate the thermal exponent yt or the droplet exponent 
9d simply by measuring the domain wall energy for sys- 
tems of various sizes instead of computing derivatives of 
the free energy or the droplet excitation energy. In fact, 
it is technically much easier to compute the domain wall 
energies than droplet excitation energies. This is one of 
the reasons why in the early stage of the study on spin 
glasses in two dimensions, many works were devoted to 
computation of the domain wall energy. We can see this 
in Tabic I. In Table I, some of previous numerical works 
for estimating one or two of the three scaling exponents, 
y t , 9s and 9o, for models with continuous, mostly Gaus- 
sian, bond distributions in two dimensions are listed in 
chronological order. 

The first estimate of the stiffness exponent was car- 
ried out through the numerical transfer matrix method 
by McMillan EP He applied the periodic or anti-periodic 
boundary condition in one direction while replicating the 
same L x L system periodically in the other direction to 
obtain a long strip, thereby mimicing the periodic bound- 
ary condition. The width of the widest strip was L = 8. 
By switching from the periodic boundary condition to 
the anti-periodic one, a domain wall is created across 
the strip. The excitation energy of such domain walls 



9 S = 0.281(5). 



(2.1) 



Bray and MooreS used the zero-temperature transfer 
matrix method in order to directly calculate the domain 
wall energy for slightly larger L's than McMillan's previ- 
ous calculation. Instead of using periodic or antiperiodic 
boundary condition for creating domain walls, they dealt 
with L X (L + 1) systems with periodic boundary condi- 
tion in the direction of L and randomly fixed boundary 
condition in the direction of L + 1. Domain walls were 
created by simultaneous inversion of all the fixed spins 
on one of two boundaries perpendicular to the direction 
of L + 1. They obtained an estimate of the stiffness ex- 



ponent onlyj-slightly smaller than Eq. (2.1). 

McMillarJlS' obtained another estimate of the stiffness 
exponent based on zero-temperature stiffness calculation 
through a heuristic optimization. An ordinary Monte 
Carlo method at a finite temperature was used for gen- 
erating many initial spin configurations. By quenching 
these spin configurations, approximate solutions were ob- 
tained. If the number of trials is large enough, i.e., if 
one generates sufficiently many approximate solutions, 
one can expect that the solution with the lowest energy 
among them coincides the true ground state with fairly 
large probability. In this way, L x L systems with pe- 
riodic boundary condition were examined up to L = 8. 
The resulting estimate turned out to be consistent with 
the above two estimates. 

Rieger et al.EJ) performed a similar calculations on 
a much larger scale with periodic boundary condition 
in both the directions. They dealt with systems up to 
L = 30 using an exact optimization method based on the 
branch- and- cut algorithm.EJ They obtained a value in 
good agreement withJEq. fl2.ip. 

Matsubara et alS3 argued that the value of the 
stiffness exponent depends on the boundary condition. 
They performed a Monte Carlo simulation on systems 
of L x (L + 1) up to L = 22. For the first system, the 
periodic and the free boundary condition was applied to 
the direction of L and L + 1, respectively. Then, for the 
second system, spins on one of the two boundaries per- 
pendicular to the direction of L + 1 are fixed as they are 
in the ground state of the first system while those on the 
other boundary are fixed opposite to the first system. 
The boundary condition in the other direction, i.e., the 
direction of L, remains to be periodic. They obtained 
an estimate of the stiffness exponent significantly larger 
than previous estimates. They argued that the stiffness 
tends to be estimated smaller when periodic boundary 
condition is imposed because it introduces an additional 
tension into a system, and this was why they used free 
boundary condition for the first system. 

In order to check if the stiffness exponent depends on 
the boundary condition, we perform a numerical calcula- 
tion of domain wall energy applying free boundary con- 
dition in every direction for the first system. For the 



Table I. Some of previous estimates of scaling exponents for the two dimensional EA model with the continuous bond distribution. 
The bond distribution is Gaussian for all the works listed here except for two items with footnotes. T m i n is the lowest temperature at 
which the computation for the largest system was performed. In the last column, TM, MC, HO, EO and FF stand for the numerical 
transfer method (TM), the Monte Carlo simulation (MC), a heuristic optimization (HO), an exact optimization (EO) and mapping 
to a free fermion problem (FF), respectively 
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(1) ... The estimate of l/^ii. (See text) 

(2) ... Periodic tiling. (See text) 

(3) ... The exponential bond distribution as well as the Gaussian distribution was used. 

(4) ... The distribution used was symmetric, consists of two continuous parts, and has a vanishing weight around J ~ 0. 

(5) ... This may not be considered as an estimate of yt- (See text) 
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Fig. 2. The domain wall energy. The typical magnitude of sta- 
tistical errors is about a quarter of the symbol size except for the 
error bar of L = 48 being slightly larger than the symbol size. 



second system we apply a fixed boundary condition in 
one direction just as Matsubara et al. did, i.e., we fix 
spins on one boundary as they are in the first system 
and those on the other boundary opposite. We carry out 
our calculation for L x L systems up to to L = 48 using an 
exact optimization method£J In Fig. 1, some randomly 
chosen examples of the domain walls are illustrated. 
In Fig. 2, the domain wall energy is plotted against the 



system size in logarithmic scale. The linearity is rather 
good and the slope is estimated as 

-6 S = 0.290(10) 

in agreement with other previous results. Here, the cited 
error has been estimated manually after doing some trial- 
and-errors in the plot. This error estimate turns out rel- 
atively conservative compared to previous estimates. In 
fact, we should be conservative in the present case be- 
cause the exponent 9s is small in absolute value which 
makes the corrections to scaling relatively large. For ex- 
ample, if we plot the excitation energy against L + a 
instead of L, where a is some constant term of order of 
unity, the estimate of 9s would change significantly. This 
sort of correction due to a constant term added to the 
system size exists, in principle, and should be taken into 
account in order to make the estimate reliable. We sup- 
pose this is one of reasons for the discrepancy between 
our estimate and Matsubara et al.'s estimate. 

We also investigate geometrical properties of domain 
walls. In Fig. 3, averaged length along the perimeter P 
and the roughness R of droplets are plotted against L. 
Here the roughness is the difference in the x-coordinatcs 
of the left-most site and the right-most site on the do- 
main wall. The obtained data for the perimeter can be 
fit well by the following scaling form 

P oc L Dlv where D w = 1.28(2). 

The estimate of D\y agrees very well with the previous 
estimate Dyy = 1.34(10)EZP for systems with the peri- 



Therefore we have 



M (T, H, L) = L d / 2 m(TL v \HL Vt+d / 2 ). 



(2.3) 



at 




Fig. 3. The averaged perimeter P and roughness R of domain 
walls. 



odic boundary condition. On the other hand, the data 
for the roughness in the logarithmic scale clearly shows 
a systematic decrease in the slope converging to unity 
around L ~ 32, which indicates a presence of a strong 
correction to scaling. However, this non-linearity can be 
easily accounted for by adding a constant to L as dis- 
cussed above. The resulting estimate is consistent with 



R oc where ( = 1 . 



2.2 



Zero-temperature finite size scaling and the thermal 
exponent y t 

As is discussed above, in the droplet argument or the 
domain- wall renormalization-group argument, we have 
only one characteristic length scale at each temperature. 
It was also concluded, based on these arguments, that all 
the exponents that characterize the scaling properties of 
various static quantities can be derived from the droplet 
exponent or the stiffness exponent. It should be noted, 
however, that we can derive this "one-parameter-scaling" 
without any special pictures or assumptions. Instead, we 
simply assume the standard finite size scaling form for 
the singular part of logarithm of the partition function, 



lo S Z sing (T,H,L) = f(TL y \HL y 



(2.2) 



where / is a scaling function and T, H, L are the temper- 
ature, the uniform magnetic field and the system size, 
respectively. From this form, we obtain for the magneti- 
zation the following expression, 

M(T, H, L) = L yh - yt m(TL yt ,HL Vh ). 

Since the model does not have a degeneracy in the ground 
state except the trivial one with respect to simultaneous 
inversion of all spins, the magnetization in the limit of 
H — > +0 should be proportional to L d / 2 . This yields a 
scaling relation 



yh = yt + 



Then, it is easy to see that the zero-temperature linear 
susceptibility at H — has the asymptotic form 



X(T = Q,H = 0,L) oc V 



(2.4) 



Similarly, the magnetization at zero-temperature has the 
form 



m(T = 0, H, L = oc) cx H 1/5 where 5 = 1 



Vt- 



There were some numerical attempts on computing 
thermodynamic quantities even before McMillan's cal- 
culation of stiffness. Such— measurements were done by 
Monte Carlo simulations^ In retrospect, their esti- 
mates were affected by insufficient equilibration in the 
low temperature region. The first direct numerical cal- 
culation of a thermodynamic quantity after researchers 
started taking into account such equilibration problems 
was carried out through an optimization method by 
Kawashima and Suzuki.© We obtained ground states 
of the systems with Gaussian bond distribution, thereby 
measuring the magnetization at a small but finite mag- 
netic field. In particular the zero-temperature linear sus- 
ceptibility for various system sizes up to L = 20 was mea- 



sured. By matching the resulting estimate with Eq. (2.4) 
the thermal exponent was estimated as 



Vt = 0.476(5). 



(2.5) 



Surprisingly, the value did not agree with the stiffness 
exponent which had been estimated rather accurately. 

From Eq. ([2.2[ ) the scaling form for the spin glass sus- 
ceptibility can be also derived. Namely, 

Xsg = 4 Tl(SiS s ) 2 ] ~ -^J^zio g z ah 



L d Xsg (TL y \HL yt 



3L d d(PH) i 

hd/2\ 



(2.6) 



where [• • •] stands for the bond configuration average 
whereas (• • •) the thermal average. The spin glass suscep- 
tibility was computed numerically through the transfer 
matrix methodic for system sizes up to L = 16 at vari- 
ous temperatures. From the best fit to the form Eq.(2.6) 
with setting H = 0, the exponent y t was estimated as 

Vt = 0.48(1), 

which agrees with the previous result obtained through 
the computation of the magnetization. 

Liangjij) also computed Xsg f° r larger systems by a 
cluster Monter Carlo method. It was at T = 0.4 that the 
simulation of the largest size (L = 128) was performed. 
The estimate of the thermal exponent was yt = 0.50(5) 
reconfirming the above-mentioned previous results. As- 
suming that the temperature dependence of the corre- 
lation length obtained by Cheung and McMillarO is 
correct down to this temperature, we can roughly esti- 
mate the correlation length at T = 0.4 to be a few tens 
of lattice constants. Therefore, Liang's calculation can 
be considered complementary to the above two previous 
calculations, in that Liang's calculation was performed 



at temperatures higher than the cross-over temperature 
whereas for the two calculations mentioned above were 
carried out at temperatures lower than that. 

The computation of the zero-temperature magnetiza- 
tion was also redone for larger systems. Rieger et al.EZP 
used an exact optimization method for computing the 
magnetization in a magnetic field for systems with peri- 
odic boundary condition. Matching the data for systems 
up to L = 50 to Eq. (2.3), they obtained an estimate of 
yt essentially identical to the previous one Eq. ([2.5]). 

Given these results on the thermal exponent together 
with the above mentioned estimates of stiffness expo- 
nent, it is now rather unlikely that the discrepancy be- 
tween them is only a numerical artifact. If they are truly 
different as the numerical estimates suggest_a simple 
domain-wall renormalization-group argument© and an 
over-simplified droplet picture. 

2.3 Droplet argument and droplet exponent 

Then, how can we elaborate the drolet picture in order 
to explain the discrepancy? A recent workS* on droplets 
may shed light on this question though the answer has 
not yet been obtained. "Droplet" is a well-known and 



Fig. 4. Droplets of 16 randomly chosen samples. 

useful concept in the study of critical phenomena. 

The droplet argument for spin glasses© based on this 
concept is one of important working hypotheses in the 
field. Droplets are collective excitations from some or- 
dered state below the transition temperature. A naive 
interpretation of the droplet theory leads to the conclu- 
sion —8s = —0d = y t . However, because of technical 
difficulty mentioned below, no direct observations of the 
droplets was carried out until quite recently. In addi- 
tion, as we have seen above, at least the simplest ver- 
sion of the argument does not work in the present case. 
Therefore, the equivalence between y t and —0d should 
be re-examined carefully. n 

Fisher and Huse defined® a droplet of scale A includ- 
ing a given site j as a cluster of spins (with i among 



them) with the smallest excitation energy that contains 
more than X d and less than (2A) d spins. The basis of the 
argument is the following scaling form that describes the 
excitation-energy distribution of droplets of scale A: 



Px(E x 



-P 



t\ o \r\ e ° 



(2.7) 



where T is some constant, #d is the droplet exponent, 
and P{x) is the scaling function which is continuous and 
non- vanishing at x — 0. Since the droplet "size", A, is 
defined to be proportional to (volume) 1 ^, it is not neces- 
sarily a spanning length of droplets because the volume 
of a droplet may in general have a non-trivial fractal di- 
mension. While it was arguedB that the droplets must 
be compact if 9d is positive, no extensive discussion had 
been devoted to the geometrical properties of droplets in 
the present two-dimensional case, where 9d is predicted 
to be negative, until quite recently. We demonstrated© 
that a typical large droplet in two dimensions occupies 
only an infinitesimal fraction of the volume of the mini- 
mal box that can contain the droplet. It follows that for 
a given droplet we have at least two essentially different 
length scales, A and the spanning length. We refer here- 
after to the fractal dimension of the droplet volume as 
D. Then, we have the relation 



D 



V = X d oc I 



where I is the spanning length of the droplet. 

Since the above-mentioned definition of droplets by 
Fisher and Huse is inconvenient from computational 
point of view, we adopted the following alternative def- 
inition. First we consider an L x L system with free 
boundary condition. A droplet of scale L is then defined 
to be the cluster of spins that has the smallest excitation 
energy among those which contain the central spin and 
contains no spins on the boundary. 

If we assume the conventional droplet argument, the 
central spin is surrounded by droplets (in the original 
definition) of various scales. Considering the fact that 
larger droplets tend to have smaller excitation energy, 
one may expect that the spanning length of the droplet 
with the new definition would be proportional to L. We 
here define 6' D by the following L dependence of the av- 
erage excitation energy of the droplets in our definition, 

[E D (L)\ kL 8 '". 

Comparing this to Fisher and Huse's definition of the 
droplet, we obtain 

— V 

because L is proportional to the spanning length of the 
droplet. 

In order to observe droplets with the new definition, 
we first computed the ground state with free boundary 
condition, and took it as the reference spin configuration. 
Then we computed the ground state with the constraint 
that the spins on the boundary are to be fixed as they 
are in the reference state while the central spin is to be 
fixed opposite, thereby forcing a cluster of spins includ- 
ing the central spin to flip. For a system with the free 



9' 

D ■ 



(2- 



boundary condition, polynomial-time optimization algo- 
rithms are available whereas for the systems with con- 
straints no such algorithm is known. In fact, the two 
dimensional spin glass problem with general constraints 
has been proven to be NP hard.Ex Therefore, we have 
employed the replica optimization!? which is a heuristic 
optimization algorithm based on the idea of renormaliza- 
tion group. j-The details of the algorithm are described 
elsewhere B cJ 

Observed droplets varies in size and shape. Some of 
them contains only one spin while spanning length of 
some others turned out to be comparable to the system 
size itself (See Fig. 4). This reflects the fact that \8d\ is 
small and the excitation energy does not very strongly 
depend on the droplet size, yielding a non-negligible 
probability for a small droplets being chosen. In the 
larger droplets, on the other hand, many handles and 
overhangs can be observed, which already suggests the 
fractal nature of the droplets. The averaged spanning 
length was found to be proportional to the system size 
as we expected. 

The averaged length, P, of the boundary of the droplet 
and the averaged volume, V were also measured. For P, 
we obtained P oc L Ds with the surface fractal dimension 



D s = 1.10(2). 
For V , we estimated the fractal dimension as 
D = 1.80(2). 



(2.9) 



The fractal dimension D of droplets is certainly smaller 
than d = 2. Thus, we concluded that the droplets at the 
critical point T = have fractal nature in the volume as 
well as in the perimeter. 

It was also found that the droplet excitation energy 
Ed(L) have a broad distribution, similar to the size 
and the shape. When rescaled with the average value 
[Eu(L)], histograms of excitation energies for various 
system sizes fit on top of each other, showing the valid- 
ity of the scaling form Eq. ( |2.7| ). In addition, we observe 
that the scaling function P(X) has a non- vanishing value 
at X = 0, satisfying a necessary condition for the droplet 
argument to be valid. 

When log Ed is plotted against log L the slope is esti- 
mated to be —0.42. From this value, together with two 
others obtained by using other measures of the droplet 
size instead of L, we obtained the estimate of the droplet 
exponent 9' D : 



9' D = 0.42(4). 



(2.10) 



Assuming Eq. (2_^8) and using the value of D in Eq. ( |2.9| ) 



(2.11) 



we obtained the following estimate of 9d, 

-6 D = 0.47(5) 

in go od agreement with previous estimates of yt such 
Eq. (|]|. 



as 



2.4 Other related estimates of scaling exponents 

As far as we know, the first important attempt to es- 
timate yt was_made through the computation of corre- 
lation length.E3 The numerical transfer matrix method 



was employed for strips of L x N where L < 11 and 
N is so large that we can take it as infinite effectively. 
The periodic boundary condition was imposed across the 
strip. Correlation lengths £|| (T, L) in the longitudinal di- 
rection, i.e., the direction of N, were measured at var- 
ious temperatures down to T = 0.15. They estimated 
the correlation length exponent as v\\ = 2.96(22). We 
should note here that this is not necessarily equal to 
v = 1/yt because of the very large aspect ratio of the 
systems. As we have seen above, a number of numeri- 
cal results suggest that domain walls and droplets of the 
same length scale correspond to different energy scales, 
and it is the droplet excitation energy, not the domain 
wall energy, that scales with the thermal exponent yt- 
On the other hand, in a long strip, it is the domain wall 
type excitations that determines the correlation between 
two points seperated by a very long distance in the lon- 
gitudinal direction. Therefore it is not very surprising 
that their estimate turned out to be close to 8s, not yt- 
Another computation based on the measurement of 
correlation lengths was performed by Huse and Mor- 
gensterr£3) for models with the exponential bond dis- 
tribution as well as those with the Gaussian distribu- 
tion. Again, the transfer matrix method was used. The 
maximum width of the strip was L = 8. Assuming the 
phcnomenological renormalization group relations, 

U ] \(T,2L) = £ ] \(T',L) 



with 



rpl 

— « 2- 1/ly n 
T 



they obtained the estimate v\\ ~ 4.2(5). 

Huse and Kot2P carried out a unique estimation of an 
exponent. They considered a system that consists of an 
infinite periodic repetition of many L x L systems identi- 
cal to each other. Since the system size is infinite, it has a 
phase transition, presumably, of the Ising type. Mapping 
the problem into a free fermion problem the eigen values 
of its row-to-row transfer matrix were obtained within 
a polynomial computational resources. Transition tem- 
peratures, then, correspond to those at which the gap 
between the largest and the second largest eigen values 
vanishes. Assuming that the typical transition temper- 
ature is proportional to L~ Vt where L is the size of the 
unit cell, they obtained an estimate yt — 0.37. However, 
again because of the existence of multiple length scales, 
it is not clear whether the exponent estimated in this 
way should be considered as yt- 

§3. Discontinuous Bond Distributions 

The critical behavior of the systems with discontinuous 
bond distributions is believed to be qualitatively different 
from the continuous distribution discussed in the last 
section because of the high degeneracy in the ground 
state. In this case the critical properties are not as fully 
understood as in the case of continuous distributions. 
What is missing is a theory corresponding to the droplet 
argument for the case of continuous distributions. 

The bond distribution in this class that is most often 




Fig. 5. The schematic phase diagram of ±J model on a square 
lattice. The thin dashed line stands of the Nishimori line, not a 
phase boundary, whereas the dotted line stands for a speculative 
phase boundary which has not yet been confirmed. 



considered is the binary bond distribution where a cou- 
pling Jij takes on the value + J and — J with probability 
p and 1 — p, respectively. For this bond distribution, the 
possibility of a finite temperature phase transition even 
in the case of the symmetric distribution (p — 1/2) has 
not been ruled outM as mentioned below. 

3.1 Intermediate Phase? 

Another problem that had been remaining to be solved 
until recently concerns the p — T phase diagram of 
the models with asymmetiicJpinary bond distributions. 
It had been suggested!^ c3'E3' that some intermediate 
phase exists in the p-T phase diagram close to the verti- 
cal part of the phase boundary (See Fig. 5) between the 
ferromegnetic phase and the other, presumably param- 
agnetic, phase. The suggestion was made based on zero- 
temperature numerical calculations. For the site-random 
model the evidence for the existence of an intermediate 
phase seemed ±a Joe even stronger than for the bond- 
random model E3' El _ 

We reinvestigatcdEZP this issue by studying the domain 
wall energy at zero temperature via the determination of 
exact ground states for large system sizes (L < 32) and 
large number of samples (iV samp i > 32768) using_a poly- 
nomial algorithm described by Barahona et alEZP We 
calculated the domain wall energy Ew for various sys- 
tem sizes with various ferromagnetic bond concentration 
p, applying the periodic or the anti-periodic boundary 
condition in one direction and the free boundary condi- 
tion in the other. The stiffness exponent 9s and another 
exponent p defined by 



[E w ] oc L" 



(3.1) 



were estimated. A positive value for p signifies the stabil- 
ity of a ferromagnetic ground state even in the presence 
of thermal fluctuations and therefore the existence of the 



ferromagnetic long range order at finite temperature.™ 
We defined p c and pi as the critical concentrations of 
ferromagnetic bonds at which [Ew] and [-E 2 ^] 1 / 2 , respec- 
tively, change from increasing to decreasing as functions 
of L. 

We hypothesized the following finite size scaling forms 
for E w and E w 



[E w ]L^ m 



MP 



(3.2) 
(3.3) 



The parameters were estimated as follows. 

= 0.896(1), = 0.77(1), = -0.19(2)(3.4) 



Pr 



p c 2) = 0.894(2), 



0.79(6), V (2) = -0.16(4)(3.5) 



The value of p^ is consistent with, most of previous esti 



(2)E2P and 0.89(10 while 

(2) 



mates such as 0.88 (2) Ja) 0.j 

inconsistent with 0.885(1).H2P The estimate of p c ' is 
larger than but marginally consistent with all the pre- 
vious estimates such as 0.86(2), H 0.850 and 0.870 jB 
while it is clearly inconsistent with 0. 854(2). c3> The co- 
incidence between and p c 2 ^ suggests the absence of 



the intermediate phase. 

It is also remarkable that not only p, 



(2) 



but also 



and tjj^ agree with the corresponding values in Eq. (3.4) 
within the statistical errors. While the agreement in p c 
already suggests the absence of the intermediate phase, 
we consider the agreement in the critical indices as an- 
other evidence for the absence of the intermediate spin- 
glass phase, since otherwise it is hardly imaginable that 
the first and the second moment of Ew show the same 
critical behavior at different values of p c . 

We also computed the stiffness for the site-random 
model. We found that pc and p^ for the site-random 
model agreed with each other and that the first and the 
second moments were scaled with the same scaling ex- 
ponents as the bond-random model mentioned above, 
though with larger corrections to scaling. 

Prior to our work, Nishimori argued© that the phase 
transition is of purely geometrical nature, which leads to 
the strictly vertical phase boundary that seperates the 
ferromagnetic phase from the paramagnetic (or glassy) 
phase in the low temperature region. If this is the case, 
the value of p c at T = must coincide that of the mul- 
ticritical point, i.e., the solid circle in the phase diagram 
Fig. 5. The most recent calculation Jor the multicriti- 
cal point was done by Ozeki and Ito,L3 who performed 
Monte Carlo simulations at finite temperatures along the 
Nishimori line to locate precisely the multicritical point 
in the phase diagram It is unsettling that previous es- 
timates of p m . c . JlSlfJ including Ozeki and Ito's estimate 
Pm.c. = 0.8872(8), are slightly smaller than our estimate 
of p c beyond the extent of errors although the difference 
is rather small. The easiest explanation is that system- 
atic errors might not be fully taken into account in the 
estimates, resulting in an underestimate of the errors. 
Another possibility is that the size dependence, Eq. (3.2) 
and Eq. (3.3), assumed in the estimation of the expo- 
nents may not correctly describe the true behavior. The 



last possibility is that the phase boundary may not be 
strictly vertical in spite of Nishimori's argument. So far, 
we have not yet identified the reason for the discrepancy. 

3.2 Relation to the Percolation Transition 

We should notice the coincidence between th^ expo- 
nent yp (or y^) in the above mentioned workEP with 
the scaling exponent y p = l/v p where v v is the correla- 
tion length exponent for the percolation tranition in two 
dimensions. Related to this finding, in a work preced.- 
ing our estimation of yjp and y p 2 \ Singh and AdleicJ 
investigated the critical behavior around the multicrit- 
ical point, and discussed its relation to the percolation 
transition. They suggested that one of the two scaling 
axes at the multicritical point is parallel to the Nishimori 
line whereas the other is parallel to the temperature axis. 
They estimated the scaling exponent which corresponds 
to our and found that the value agrees with the per- 
colation exponent. Since, according to Nishimori's argu- 
ment, the temperature plays no essential role along the 
vertical part of the phase boundary, we could argue that 
the phase transitions across the boundary are character- 
ized by the same scaling exponent independent of the 
temperature, although we are not sure if the argument 
should be valid at the two end points as well. 

3.3 Finite- Temperature Phase Transition at p = 1/2? 
In the case of the bond-random ±J model with p = 

1/2, the majority of researchers tend to think that there 
is no phase transition at a finite temperature. This 
generalj-belifif is based on results of Monte Carlp_f;imu- 
lations,ErE3 high-temperature series expaaswnjp' and 
the estimates of the domain wall energy .BUSES The 
data from Monte Carlo simulations are not available at 
very low temperature, e.g., below T = 0.4 in Bhatt and 
Young's simulation,© which makes it difficult to exclude 
the possibility of the transition at a temperature smaller 
than 0.4. Although the results of high-temperature series 
expansion also suggested the absense of the finite temper- 
ature phase transition, it was not very conclusive. As for 
the calculations of the domain wall energy, the "stiffness" 
exponent turned out to be even smaller in absolute value 
than that for the continuous distributions. Therefore, it 
was not very clear if 9s is really negative. For example, 
the data in Cieplak and Banavar's papertS clearly shows 
a systematic positive curvature in the plot of the domain 
wall energy versus the system size. In short, so far we 
have not obtained very convincing evidences supporting 
the general belief, i.e., the absense of the phase transition 
at a finite temperature in two dimensions. 

On the other hand, results of Monte Carlo simulation 
at lower temperatures were reportcdES indicating a tran- 
sition at T ~ 0.24 for the ±J model with p = LZ2. A 
finite temperature phase transition was suggesteclEf also 
for the model with another discrete distribution where a 
bond variable takes on J or — aJ_where < a < 1. In 
addition, 0$ — was suggested^ based on estimates 
of the domain-wall energy. These results are consistent 
with a finite temperature phase transition for which the 
low-temperature phase is only marginally or weakly or- 



dered, i.e., the two-point spin-spin correlation function 
decreases algebraically as ^Junction of the distance. 

On this issue, we foundEZP that the stiffness decreases 
systematically but it does so extremely slowly in the case 
of p — 1/2. The obtained data was consistent with the al- 
gebraic decrease as a function of the system size with the 
stiffness exponent 9s = —0.056(6). However, we should 
note that the quoted error here only includes the esti- 
mated statistical error that is obtained in the standard 
procedure of the method of least squares. Because of the 
very weak dependency on the system size, the systematic 
error is more important here than the statistical error. 
Since the systematic error is difficult to estimate, iithcr 
scenarios, such as convergence to a finite constant© or 
logarithmic decay, are equally probable. If cither one 
of these scenarios is correct, we cannot determine based 
solely on calculations of the stiffness exponent whether 
the long-range order persists at a low but finite tempera- 
ture. At the moment, what we can conclude is only that 
the low-temperature phase is only weakly ordered even if 
the phase transition takes place at a finite temperature. 

Very recently, Kitatani and SinadaEf performed an 
interesting numerical calculation of spin glass suscepti- 
bility and its Binder ratio. Based on the invariant prop- 
erties under the gauge transformation, they constructed 
a new method for computing cumulants of Edwards- 
Anderson order-parameter at any temperature using a 
novel Monte Carlo method in the bond configuration 
space rather than in the spin configuration space. Their 
results of the finite size scaling analysis, appeared to sug- 
gest the existence of a phase transition around T = 0.3. 
However, based on the data analysis with corrections to 
scaling, they concluded that the scenario of the zero- 
temperature phase transition could explain their data as 
well. 

§4. Concluding Remarks 

We have reviewed recent developments in study of spin 
glass models in two dimensions. It has become clear in 
the case of continuous bond distributions that the critical 
behavior is dominated by fractal droplets in contrast to 
the compact droplets in conventional droplet argument 
for higher dimensions. 

The disagreement between the droplet exponent and 
the stiffness exponent remains unsettling. At least, how- 
ever, the present result suggests that a domain wall of 
a L x L system may not necessarily be considered as 
an object of scale L because of the existence of different 
ways for measuring the size of the same object. In other 
words, at present we do not know which way of mea- 
suring we should use for comparing the size of a domain 
wall with that of a droplet. 

The models with discrete distributions are relatively 
poorly understood. It is an important and challenging 
future problem to determine whether a finite tempera- 
ture phase transition takes place or not by larger scale 
numerical calculations. 

Although it is not yet clear whether fractal objects 
similar to the present two-dimensional case can be de- 
fined at the critical point in higher dimensions, it is cer- 
tainly worth thinking about such a possibility. The es- 



sential difference wouid be that rj is not in general zero 
in higher dimensions whereas 77 = for continuous dis- 
tributions making the theory simpler in the latter case. 
Provided that the ± J model in two dimensions is crit- 
ical at T — 0, the same remark may apply to this dis- 
crete model since rj is non-vanishing there, too. There- 
fore, we may obtain some new perspective for both the 
finite-temperature transition in three dimensions and the 
zero-temperature critical behavior for discrete distribu- 
tions in two dimensions by generalizing the concept of 
zero-temperature fractal droplets to the cases with non- 
vanishing 77. 
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